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Abstract 


A 3D micro-scale model is developed to simulate the transport and electrochemical reaction in a composite cathode. This model takes into 
account the details of the specific cathode microstructure such as random pore structure, active TPB (three phase boundary) site distribution, 
particle size and composition and their interrelationship to the charge transfer and mass transport processes. Especially, the pore structure and 
mass diffusion were incorporated into this model. Influence of the microsturcture parameters on the performance was investigated by numerical 
simulations. 

Simulation shows that the cathode porosity should be in the range of 0.25-0.45 for the optimized performance. A larger thickness is in favor of 
increasing the effective reaction sites and reducing total specific resistance. However, the thickness needs to be confined below a certain thickness 
in order to prevent lower utilization and excessive concentration loss. The fine particles contribute to increase the TPB length and thus decrease the 
activation overpotential. On the other hand, the extremely small particles might lead to an excessive diffusion resistance. The model was validated 
based on the continuum model and experimental data to determine its accuracy in predicting the cathode performance. The results generally agree 
with the experimental data from the literature. However, as compared with the continuum model, the predicted total resistance is slightly higher 


for a thin electrode. 
© 2007 Elsevier B.V. All rights reserved. 
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1. Introduction 


Solid oxide fuel cells can provide clean and highly efficient 
energy utilization in wide applications from small stationary 
power units to larger scale capacities. Many investigations on 
the loss mechanisms by impedance spectroscopy as well as the- 
oretical analysis have revealed that the cathode governs the main 
part of the losses and thus the improvement in electrochemical 
behavior strongly depends on selected material and microstruc- 
ture of the cathode, etc. In the case of pure electronic-conducting 
electrode materials like metals or some perovskite type oxides 
(LSM), electrochemical reactions are almost confined to the 
triple phase boundaries (TPB, gas/electrolyte/electrode) [1,2]. 
The transport of oxide ions within the electrode materials 
is advantageous concerning the number of possible reaction 


* Corresponding author. Tel.: +1 352 392 9607; fax: +1 352 392 1071. 
E-mail address: jnchung @ufl.edu (J.N. Chung). 


0378-7753/$ — see front matter © 2007 Elsevier B.V. All rights reserved. 
doi: 10.1016/j.jpowsour.2006. 12.058 


pathways. Therefore, high performance electrodes should be 
composed of either a composite consisting of an electronic 
and an ionic conducting phase or a mixed conducting metal 
oxide to expand the active reaction area into the electrode vol- 
ume. The performance of a composite electrode is determined 
by the following factors: (1) conductivities of pure ionic and 
electronic conductors, (2) specific surface area available for elec- 
trochemical reaction (or the size of TPB), (3) microstructure 
of electrode, i.e., the compositional volume ratio, porosity and 
particle size distribution, and (4) path of gas diffusion. For exam- 
ple, Kim et al. [3] indicated that polarization resistance of the 
LSM-YSZ composite electrodes is closely related to the con- 
nectivity of LSM—YSZ particles and depends on their size ratio. 
The composite electrode showed low electrochemical activity 
without YSZ connectivity. However, the polarization resistance 
decreased abruptly when the YSZ connectivity is present. Thus, 
optimization of conductivity and polarization resistance of com- 
posite electrodes for a SOFC in terms of the microstructure will 
undoubtedly gain substantial understanding. 
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Nomenclature 

BS diffusion matrix 

Č mass concentration vector (kg m~?) 

De diameter (m) 

DP, binary mass diffusivity (m? s7!) 

DS effective Knudsen diffusion coefficient for com- 
ponent i (m? s7!) 

Fy Forchheimer coefficient 

itr exchange current density of cathode (Am7! or 
Am?) 

I current flow (A) 

j diffusion mass flux vector 

Ly, Ly, Lz dimensions of electrode (m) 


Mi molar weight of component i (kg mol!) 

Ny, Ny, Nz size of lattice in x, y and z directions 

p pressure (Pa) 

Po, partial pressure of oxygen at reaction site (atm) 
Tp aveage pore radius (m) 

R rate of consumed mass (kg m7? s7!) 

Rg universal gas constant, 8.3143 J mol~! K7! 

Riot total resistance (Q) 

T temperature (K) 

Xi molar fraction of component i 


Greek letters 


B transfer coefficient 

ô neck of diameter (m) 

E porosity 

(o) potential (V) 

y contact angle (°) 

Nact activation overpotential (V) 

Neon concentration overpotential (V) 

Kio, Kel ionic conductivity, electronic conductivity 
(27! m7!) 

u viscosity (m? s7!) 

oO electrical or ionic conductance (Q7!) 

wij characteristic length 

Qp collision integral 

E coordinate of the pore 

Aé;_; distance for any two neighbor particles (m) 

£ neck circumference (m) 

Subscripts 

el electronic 

io ionic 


However, the microstructure of an electrode is not regular 
and thus all mass and energy transport within electrodes take 
place in a disordered, unorganized medium. This random char- 
acteristic makes modeling electrode become extremely difficult, 
because the accuracy of a model is dependent on how to describe 
this random geometry. One of the representative publications 
is Costamagna et al.’s continuum model [4] based on the per- 
colation theory and coordination number theory [5—10]. This 


statistical method gives a good starting point for composite elec- 
trode simulation due to simplicity and explicity; however, some 
inherent assumptions and oversimplifications limit the appli- 
cation range of this type of continuum-model: (1) the mass 
resistance in a porous medium is completely neglected from the 
assumptions of faster diffusion and smaller particle size. How- 
ever when the thicknesses of electrode or current density are 
greater, the overpotential due to mass resistance is not negligi- 
ble. (2) Since only sample-spanning electrolyte and electrode 
cluster are accounted for in this model, the theory is not valid 
outside the composition regime in which isolated clusters are 
present and has lower ratio of thickness to particle radius. (3) 
The mechanism of electrochemical reaction was not investi- 
gated. In order to further apply the continuum-model to actual 
SOFC electrodes, some modifications are proposed by Chan et 
al. [11]. The primary improvement is the electrochemical reac- 
tion in which the exchange current density is not constant as in 
previous work. Jeon et al. [12] extended the original model to 
an actual multilayer electrode. 

Obviously, the most natural modeling methodology is 
directly based on the actual microstructure of such composite 
electrodes. In principle, any disordered system can be mapped 
onto a random network of bonds connected to each other at sites 
or nodes of the network. The network models incorporate pore- 
scale descriptions of the medium and the physics of pore-scale 
events [13-20]. This method represents the type of approach that 
best predicts the most important experimental features of com- 
posite electrodes. However, their models suffered from some 
common weaknesses. These models overemphasize the effect of 
selected parameters on the performance and do not comprehen- 
sively evaluate the overall performance. Further, the influence 
of pore distribution, thermal factor and current distribution are 
rarely considered. 

Therefore, reliable models capable of handling realistic 
microstructural features and mass/heat transport are essen- 
tial for establishing a quantitative relationship between the 
microstructure and the electrical properties in composite elec- 
trodes. Accordingly, a micro-scale model considering more 
detailed mass and energy transport processes will be developed 
in the present paper. A novel 3D random walk model will be 
applied to optimize the performance in a porous composite cath- 
ode resulting in a complete description of the electrode structure 
as well as the processes occurring therein. The pore distribution, 
composition, mass diffusion and current density distribution, etc. 
will be completely integrated into the micro-model. 


2. Physical description and assumptions for 
micro-model 


The networks of porous cathode are complex, and in prac- 
tice it is necessary to adopt simplified models which reproduce 
the essential features, but ignore many of the details, of a real 
network. In the present model, an idealized 3D microstructure 
is constructed as shown in Fig. | (only 2D view is given). The 
microstructure takes the form of a face-centered cubic lattice 
(FCC). This kind of structure represents a dense packing of par- 
ticles of the same size (i.e., electrolyte, electrode and pores). A 
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Fig. 1. 2D view of random network model. Red rectangle denotes the active TPB; 
black circle, electronic conductor; blue circle, ionic conductor; white empty, 
pore (the circles have been removed for clarify). Three paths along the three 
phases for gas (black), O7 (red) and e~ (white) are represented by curves. (For 
interpretation of the references to color in this figure legend, the reader is referred 
to the web version of the article.) 


practical cathode has a random topology, i.e., its local coordi- 
nation number, which is the number of bonds connected to a 
site, is a chaotic and almost random variable. Sunde [18,19], 
Zalc et al. [15], Constantinides et al. [21] and Ciobanu et al. 
[22] employed the arbitrary resistor network with non-uniform 
grain size and shape to represent the porous structures. How- 
ever, some studies [23,24] have demonstrated that as long as 
the average coordination number of the topologically disordered 
network is equal to the coordination number of a topologically 
regular network, transport and main properties are identical for 
the two network systems. Therefore, a 3D regular lattice is used 
here. The spatial distributions of grain particles and pores are 
still random. In Fig. 1, blue, black, and white sites represent 
the electrolyte, electrode and pore particles, respectively. The 
electrons from inter-connector through electrolyte clusters and 
the oxygen through contiguous pore passage are supplied to the 
active TPB (explained later), where the oxygen ions are gen- 
erated by electrochemical reactions. Then the oxygen ions are 
transferred to the electrolyte membrane by electrode clusters. In 
developing this model, several assumptions and definitions are 
made and summarized as follows: 


(1) Not all of pores are “open” for reaction. Some of pores may 
be inaccessible to one or more of the species so that the 
effective TPB length and effective diffusivity are reduced. 
When the porosity is small, the pore forms small and isolated 
clusters. As the porosity is increased enough, there would 
be some clusters spanning the network from one end to the 
other. This kind of pore clusters is termed as “transport pore” 
cluster, since it allows the gas to transport across the entire 
medium. A pore belonging to a “transport pore” cluster is 
called a “transport pore”. The second kind of pores, which 
is called “dead-end pore”, is interconnected only from one 
side (current collector). Although these pores can often be 
penetrated, their contribution to transport is negligible to 
transport. The third kind of pores is “isolated pore”, which 
is totally isolated from its neighbors of “transport pore” 


and cannot contribute to transport of mass cross the porous 
medium. Here, this common realistic situation is consid- 
ered and pores are randomly distributed and identified as 
“transport”, “dead-end” or “isolated” pores. 

(2) Three types of clusters in a solid phase, i.e., “percolating 
cluster” (I), “singly connected cluster” (II) and “isolated 
cluster” (III) are identified. The first kind of clusters, i.e., a 
percolating cluster, is that the same type of particles (ionic or 
electronic particles) touches one another to extend through 
the entire electrode. Under such a condition, a good con- 
ductivity is reached. The second is that a cluster is only 
connected to its corresponding bulk phase type (i.e., inter- 
connector or electrolyte). The last kind of clusters, “isolated 
cluster’, is completely isolated from its corresponding bulk 
phase. The existence of this kind of clusters increases the 
polarization resistance in the electrode. The distribution of 
these clusters in electrodes is strongly influenced by the 
particle ratio, the volume ratio and sintering process, etc. 

(3) Gas is considered to be an ideal two-component mixture of 
O7 and a dilute gas (N2). Within the pores, there are three 
distinct diffusion mechanisms for mass transfer [25]: molec- 
ular diffusion, Knudsen diffusion and surface diffusion of 
adsorbed molecular species along the pore wall surface. The 
molecular and Knudsen diffusion are dominant mechanisms 
at the cathode and will be taken into account in this model. 

(4) Only the dc part of the effective electronic conductivity is 
considered. The ac part of the conductivity due to double 
layer effect (i.e., capacitance effect) is neglected. 

(5) The cathodic reaction site is limited to the vicinity of an 
“active bond”. This concept is firstly introduced by Alder 
et al. [26]. The bond between two neighboring sites, one 
of which is an electrolyte and the other an electrode site, is 
active if both of the sites belong to a “percolating cluster” or 
a “singly connected cluster” and also the bond has a neigh- 
boring empty site (pore), which belongs to the “transport 
pore” cluster or “dead-end pore” cluster. 

(6) Isothermal and steady-state operation. 


3. Derivation of numerical model 


According to the conservation principle, the current flux 
entering a control volume should be equal to those leaving this 
volume plus the stored current within it. Therefore, the govern- 
ing equations for charge transport in the solid phase (two phases) 
and gas species transport in the void space, respectively, are: 


V - (cep Vd) = 0 (1) 
V- Jot + R= 0 (2) 


Here ¢ and Jiet are the electrical potential within solid par- 
ticles and total mass flux vector, respectively. R is the rate of 
consumed mass due to an electrochemical reaction. Ceff is the 
effective conductance. The values of conductance in the solid 
phase are finite and zero in the gas phase, respectively. Actually, 
Eq. (1) is applied to both of the two solid phases and Eq. (2) is 
valid for pore space. The two equations are coupled by electro- 
chemical reaction on the active three phase boundary. Hence, the 
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Fig. 2. Schematic of electrical network. 


idealized random microstructure can be considered to split into 
two sub-networks: an electrical network for the solid phase and 
a pore network on the void space. The species concentrations 
and chemical reaction rate (Eq. (2)) will be calculated within 
the pore network, whereas the charge transfer, resistances and 
overpotentials (Eq. (1)) will be analyzed within the electrical 
networks. 


3.1. Electrical network 


Fig. 2 illustrates a two-dimensional electrical network. Nodes 
(i, j, k) represent solid particle sites. Between each pair of solid 
particles, the current is conducted through a resistance with the 
“effective conductance Oe”. There are three kinds of conduc- 
tances: the conductance between the electrolyte particles Cio-io, 
the conductance between the electrode particles o¢.; and the 
conductance between two different particles oio-e1. Here, sub- 
script “io” refers to the electrolyte and “el” to the electrode 
particles. These conductances depend on the neck formed by 
the contact between two particles. A neck of diameter 26 is 
formed between two contacting spheres with the same radius r 
and thus the neck circumference, £, is equal to 2776. The conduc- 
tance between the same types of particles can be approximated 
as 


Kiolio—io Keltel—el 


Oio—io = ae `’ Oel—-el = =a (3) 
where Kio and Ke are the ionic conductivity and electronic con- 
ductivity, respectively, which are dependent on the temperature. 
The conductance Ojo] between different types of particles is 
given below based on linear electrochemical kinetics [18,19]: 


: + : + Zi (4) 
o-a = 
os 20¢l-el 20i0-io Op 


where øp is the polarization conductance. For an inactive TPB, 
Op is set as infinity. While for an active TPB, op depends on the 
geometry of a TPB and the activity of a cathodic reaction. In 
previous works [12,18], op was assumed to be a constant and 
proportional to the polarization conductivity «p and three phase 
boundary length rpg. However, Sunde [18] has implied that his 
model has significant error due to neglecting the effect of gas 
diffusion on the electrochemical reaction and kp determined by 
a simplified geometry, not the electrolyte—electrode contacts in 
the composite. Thereafter, a new approach is used in this paper 
to evaluate this value considering the effect of gas diffusion, 


temperature, chemical reaction and TPB geometry. nact is the 
activation overpotential and expressed by 


Nact = Pel — Pio — Adconc (5) 


So the polarization conductance corrected for ohmic drop can 
be expressed by 


Nact T Pel — Pio — Adgconc =. 
o= |- = ; (6) 
iL TPB iL TPB 


here £7pp is still equal to 27rd and also can be given by 


erp = 28 = 27r sin” (7) 


where r is the radius of the particle and y is the contact angle 
between ionic and electronic conductors. Chan et al. [11] found 
that this angle can be reasonably estimated as 60°. The concen- 
tration overpotential is therefore calculated as [27] 


RoT (re 
A@cone = a In pi (8) 
02 


where Po, and Po, are the partial pressures of oxygen at cath- 
ode/channel interface and reaction site (TPB), respectively. 

In Eq. (6), there is still one parameter, the exchange transfer 
current, itr, left to be determined. As a matter of fact, oxygen 
reduction reaction has been recognized to involve a series of 
consecutive elemental steps including, i.e., dissociation, charge 
transfer, surface diffusion of oxygen intermediate species, and 
incorporation of oxygen ions into the YSZ electrolyte. The 
current density depends on the partial oxygen pressure, rate 
constants in elementary reaction steps, detailed experimental 
condition and microstructure of electrodes, etc. It is difficult to 
evaluate this value theoretically; instead, fitting from experimen- 
tal data is generally used [11]: 


, 2.05 x 104 
i = f(T, poy, act) = 7874 x Pa. x exp <a 
x {exp( fMact) — exp(— fnact)} (9) 


here po,is the partial oxygen pressure, and T is the absolute 
temperature. f= F/RgT, where Rg is the gas constant and F is the 
Faraday’s constant. Finally, the current conservation equation 
(Eq. (1)) can be replaced by the following discretized form: 


Soh) = >) aij(¢i — 6)) = 0 (10) 
J j 


where ¢; and @; are the potentials for particle i and its neighbors 
j, respectively. I;>; denotes the current flowing from site 7 to its 
neighboring site j. oj; is the conductance between sites i and j. 
By iteratively solving Eqs. (3)—(9), the distribution of potential 
in the composite cathode can be obtained. Once the potential 
distribution is determined, the total resistance for the cathode is 
estimated as 


Prop — Poot 


Riot = 
I total 


(1) 
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It should be noted that the total resistance includes all three 
overpotentials, i.e., ohmic loss, activation loss and concentration 
loss. If the polarization conductivity is set to be infinite (op = 00), 
the total ohmic resistance Ronm will be easily evaluated. Then, 
the total polarization resistance Rp is obtained by subtracting 
Rohm from Riot, 1.€., Rp = Rtot — Rohm. 


3.2. Diffusion and reaction rate 


There are three diffusion mechanisms that exist within the 
voids. The total flux of any species is obtained by combin- 
ing the separate contributions. Because the bulk and Knudsen 
diffusion are dominant mechanisms in the composite cath- 
ode which has larger particle sizes (0.1-5 wm), therefore the 
surface diffusion is neglected in this model. The total mass 
transfer in a single pore is a sum of diffusive and viscous flux 
(kg m~? s7! ): 


Jot = Jaitt + Jis (12) 


The first term on the right hand side (Jait) is the mass flux 
due to bulk and Knudsen diffusion, which can be expressed by 
a N-dimensional matrix notation: 

3 -> -1 d 

Jaitt = —[B° fp, O] — (13) 
dé 

where C and € are the mass concentration vector and coordinate 

of the pore, respectively. The N — 1 dimensional square matrix 

[B°] has the elements: 


N-1 
Xi 1 Xk 
Busey =--> By=-—=t+ — 14 
GAD D}; ii DK a, DP, ( ) 


here X; represents the mole fraction of component i, Di, the 
binary diffusivity of components i and j, and DE is the effec- 
tive Knudsen diffusivity of component i. N is the total number 
of components. Both of the two diffusivities can be calculated 
as 


2 [8R;T 
DES a apn (15) 
L 
T3(M; + M)/M:M 1” 
D? = 1.86 x 10-7! ie MM (16) 


po, 2D 


where æw;j is the characteristic length, and Qp is the collision 
integral. The viscous flux in Eq. (11) is given by the following 
expression: 


: s i dC 
Jag C E Rp eo 17 
VIS 8u dé ( ) 


here u is the viscosity and T is the temperature. Ctot is the total 
molar concentration and given by 
N 
Ci 
Co = D> (18) 


1? 


The volumetric source term R (kg m7? s7!) due to electro- 
chemical reaction is related to the transfer current density itr 
(Eq. (5)), TPB length and pore size, which can be expressed 
by 


it lTPB 


E= ayn AF ” 


It is apparent that only the gas phase site having active 
TPB has the source term R. The pore structure can be treated 
as a three-dimensional pore network, and the single pores 
are connected at the nodes of the network. The flux of each 
species that enters a node must be equal to the flux leaving 
the node plus the accumulative mass. Therefore, at the inner 
nodes of the network, an equation similar to the Kirchhoff’s law 
holds: 


SUV Jinn = Ri (20) 
ie¢M 

where M denotes all the pores that are connected to node Z, and 
— Ji; are fluxes leaving or entering node /. Actually, Eq. (20) 
is the discritized form of Eq. (2). For a two binary component 
system (O2 and N2), Eq. (20) can be rewritten in the following 
expression by combining Eqs. (13) and (17): 


y o 2g et a De 
ea Se Aii (16/3)73 


S 1 T 1— &«Xo,,1 , Coa, i E Co,,1 
an Ds, Dos,No Agj-1 
D: R 21) 
-ar E a 
16/33 
where 
M. 1/2 
a=1— (Fe) (22) 
N2 


A£; is the distance between any site J and any neighboring 
pore i. Xo,,7 is the molar fraction of oxygen at site /. In this 
model, each pair of any two neighboring sites are assumed to 
be connected by a cylindrical channel with length Aé;_; and 
diameter De. Finally, Eqs. (10) and (20) have to be solved 
simultaneously by a finite difference scheme for the entire 
network. 


4. Reconstruction of microstructure and boundary 
conditions 


Fig. 3 shows a typical constructed microstructure for the 
24 x 24 x 18 cubic FCC lattice. The size of the lattice is 
represented by Ny, Ny and N; layers in the x, y, and z direc- 
tions, respectively. The z-direction is taken as the primary 
direction of current flow, and the dimension in this direc- 
tion, is by definition equal to the thickness of electrode, L,. 
The radius of a particle is in the range from 0.1 to 5pm. 
Once the microstructure is constructed, the connectivity and 
the different clusters need to be identified. For this purpose, 
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Fig. 3. Reconstruct of 3D random microstructure. Model parameters: ¢=0.4, 
et:Pio = 1:1, 24 x 24 x 18. Blue (0), isolated particle; cyan (10), transport pore; 
green (20), transport electronic particle; yellow (30), transport ionic particle. 
(For interpretation of the references to color in this figure legend, the reader is 
referred to the web version of the article.) 


the Hoshen—Kopelman algorithm [16,28] is used to identify 
the “percolating cluster”, “isolated cluster”, “dead-end clus- 
ter” and “active TPB”. All of the quantitative results will be 
obtained as an average over six different random microstructures 
with different seed numbers. Fig. 4(a) gives the local percent- 
age variation for three different sites. These data suggest that 
the pore, electronic and ionic sites are randomly distributed 
across the thickness and the algorithm used to describe the 
microstructure is unbiased. All of the electrochemical variables 
in the z-direction are of interest, so the edge effect of x- and 
y-dimension on the electrochemical behavior should be mini- 
mized. 

Fig. 4(b) shows the dependence of the effective porosity on 
the x and y dimensions. It shows that when the total number of 
cells in the x and y directions is increased from 24 to 36, the 
transport pore volume fraction varies only slightly. Therefore, 
the number of lattices is chosen at 24 that is considered to be suf- 
ficient in the x and y directions, so itis employed for the following 
calculations. Also, it is noted that when the natural porosity is 
larger than 0.3, the effective porosity is almost identical to the 
natural porosity. However, when the natural porosity is larger 
than 0.5, the volume fraction of “transport” solid particles (elec- 
tronic and ionic) dramatically decreases since the “isolated” and 
“dead-end” solid clusters dominate the whole electrodes. There- 
fore, it is concluded that 0.3 is the lower bound of the porosity. 
In practical applications, the porosity is always in the range of 
0.35-0.45. 

The boundary before the first z layer represents the equipo- 
tential surface of the current collector (metal), which means that 
all bonds connected with it are linked to the metal and gas trans- 
port channels and the potential is set as ġo. The boundary after 
the last z layer contacts the solid bulk electrolyte membrane, 
which stands for another equipotenal surface and is set to be at 
zero voltage. In the x and y directions, symmetric boundary con- 
ditions are employed. The boundary conditions are summarized 
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Fig. 4. Statistics of random microstructure (model parameters: ¢=0.4, 
det:Pio = 1:1): (a) local variation in three kind of sites along the thickness; (b) 
effects of cross-sectional mesh resolution on the transport pore volume fraction 
and transport solid volume fraction. 


as below: 
At x=0, y=0, x= Ly, 
dCo ao) 
=Ly:—2=0, —=0 23 
7 7 on on (23a) 


At the electrode/metal interface 


z=0: Co, = Co, ea = po, I= hot (23b) 
At the electrolyte/cathode interface 
aC 
z=Lz;: Z =0, dio =0 (23c) 
on 


Note that in Eq. (23b), only the potentials for electrode par- 
ticles are set as finite values, and the potentials for electrolyte 
particles should be considered as “inner nodes” and obtained by 
solving the governing equations. Similar cases hold for boundary 
condition (23c). The input parameters, including the transport 
and kinetic parameters are extracted from literatures [4,11,27] 
and listed in Table 1. The conductivities for LSM and YSZ are 
considered to be dependent on the temperature. Air (79% N2 
and 21% Oz) is delivered to cathode. 
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Table 1 

Model input parameters for the baseline case 

Parameters Value Parameters Value 

Air inlet temperature (K) 1100 Contact angle, y 60° 

LSM conductivity (S m7!) (8.855 x 10’/T) x exp(—1092.5/T) YSZ conductivity (S m7!) 3.34 x 10° x exp(— 10300/T) 
External voltage difference (V) 0.2 Bulk oxygen concentration (kg m7?) 0.0825 

Ly x Ly x L 24 x 24 x 18 Particle diameter (wm) 2x 107° 

Porosity 0.4 


5. Numerical results and discussions 
5.1. Contours of potential and oxygen concentration 


To study the potential distribution, the mass diffusion and the 
significance of overpotentials around the overall cathode region, 
we first select a microstructure with the seed number “2” and 
“3” as the baseline case. Although results obtained from this 
random microstructure is different from that based on an aver- 
age of several microstructures with different seed numbers, the 
errors are so small that will not significantly change the values 
of potential and oxygen concentration. In this case, the exter- 
nal voltage difference is set as 0.2 V, the porosity is 0.4 and the 
composition ratio of electrode phase is 0.5. The average cur- 
rent density is calculated as 0.806 Acm~*. The 3D contours 
of potential are plotted in Fig. 5. It turns out that the poten- 
tial in the ionic phase gradually decreases from about 0.11 V to 
zero across the electrode thickness (z-direction). For electronic 
phase (Fig. 5(b)), however, the variation in potential is almost 
negligible. These distributions attribute to the higher conductiv- 
ity for LSM and very lower conductivity for YSZ. Such larger 
potential difference between LSM and YSZ is just caused by 
the electrochemical reaction in the vicinity of TPB. Distribu- 
tions of corresponding average cross-sectional exchange density 
and polarization conductance through the cathode, normalized 
to the current and the polarization conductance at the first layer, 
are displayed in Fig. 6. In the present study, the exchange current 
is modeled as a function of activation overpotential, and other 
operation parameters, which is different from the assumption 
of constant value in Sun’s study [18,19]. The exchange current 
density monastically increases and reaches the peak value at the 
bottom of electrode. This tendency indicates that the strongest 
electrochemical reaction takes place on the region near the bot- 
tom end of the cathode, where the larger potential difference and 
maximal number of active TPB exist. For pure phase cathode, 
the electrochemical reaction is only confined to a very narrow 
zone, i.e., the electrolyte/cathode interface; while for compos- 
ite cathode, the ionic phase extends to the bulk LSM phase and 
hence increases the reaction sites. Also, we found that the higher 
active TPB number, the larger the polarization conductance is. 

Fig. 7 compares the oxygen concentration under lower 
and higher potential differences. Under the lower potential 
difference, the oxygen concentration is fairly uniform, decreas- 
ing from 0.08259 kg m~? to roughly 0.075 kg m~? across the 
electrode thickness. While for the case of a higher potential dif- 
ference, i.e., 0.4 V, the oxygen gradually decreases to a lower 
level of oxygen concentration due to a larger oxygen con- 


sumption rate. The uniform oxygen distribution on a x-y plane 
validates that the percolating cluster dominates the whole region 
and makes the oxygen well supplied to the active TPB. 


5.2. Activation and concentration overpotentials 


Fig. 8(a) illustrates the profiles of average ohmic and acti- 
vation overpotentials with respect to different external voltages. 
An increase of voltage difference (i.e., average current density), 
the sum of ohmic and activation loss also increases. Although 


(a) 


Fig. 5. 3D contours of potential (unit: V): (a) for electronic particles; (b) for 
ionic particles. 
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Fig. 6. Profiles of local relative exchange current density and relative polariza- 
tion conductance. 
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Fig. 7. 3D contours of oxygen concentration: (a) lower potential difference; (b) 
higher potential difference. 
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Fig. 8. Profiles of overpotential distributions: (a) ohmic and activation loss; (b) 
concentration loss. 


both of the numbers of active TPB on the gas/cathode and 
electrolyte/cathode are close to the maximum values, the over- 
potenial next to the electrolyte/cathode surface, i.e., Z/Lz = 1, is 
still much higher because of higher polarization potential differ- 
ence in this region. Actually, the activation loss contributes to the 
major portion of the overpotential. The corresponding average 
concentration loss is given in Fig. 8(b). It indicates that this loss 
is negligible compared with the activation loss under the lower 
current density. However, when the current density is enhanced 
due to a larger potential difference, the concentration loss can- 
not be neglected any more. For example, when the potential is 
0.4 V, the maximal concentration loss is about 0.0325 V that is 
responsible for 9% of the total loss. 


5.3. Optimization of cathode design 


5.3.1. Porosity effect 

Fig. 2 implies that the cathode should have the best perfor- 
mance with the porosity exceeding 0.3 in terms of the effective 
porosity. But how the porosity exactly influences the mass 
transfer and electrochemical behavior still needs to be further 
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Fig. 9. Effect of porosity on performance: (a) average local oxygen concentra- 
tion; (b) average reaction rate. 


investigated. The average oxygen concentration distribution and 
reaction rate as a function of porosity are shown in Fig. 9(a) and 
(b), respectively. The porosity is varied from 0.1 to 0.6 while 
keeping other parameters the same as that listed in Table 1. 
Clearly, when the porosity is at a relatively lower level, i.e., 
0.1 or 0.2, the oxygen penetration depth is considerably small, 
only 25% of the cathode thickness for £ = 0.1, which will lead to 
a larger concentration overpotenital. The peak value of reaction 
rate occurs near the gas/electrode interface. This is due to the 
fact that an “isolated pore” dominates the whole electrode and 
blocks the gas supply to the TPB, and then substantially reduces 
the gas penetration depth. Therefore, the availability of oxygen 
supply and reaction sites is confined to a very narrow zone next 
to the gas/electrode interface. Nevertheless, when the porosity 
is equal to 0.3, the reaction rate begins to have a continuous 
tendency to increase along the thickness and also possesses the 
highest reaction rate. When the porosity is enhanced to 0.4 or 
more, the cathode performance deteriorates due to insufficient 
reaction sites. 

Fig. 10(a) shows the statistical characters of the active TPB 
with respect to the porosity. It is observed that when the porosity 
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Fig. 10. Profiles of number of TPB and specific total resistance with respect to 
the porosity: (a) total number of active TPB; (b) specific total resistance. 


ranges from 0.25 to 0.45, the percentage of active TPB number 
lies between 77% and 91%, which means the most of transport 
ionic or electronic particles are available to the electrochemi- 
cal reaction. Therefore, in light of the above results, it can be 
concluded that the concentration loss can be neglected and the 
cathode performance reaches the best when the porosity is in 
the range of 0.25-0.45. In Sunde’s study [18], it is assumed that 
all of solid particles are always available for reaction and the 
oxygen supply is not limited in the vicinity of a TPB. There- 
fore, their results are only valid for the porosity range in the 
middle level and also in general the resistance is underestimated 
in their study. The corresponding total resistance excluding the 
mass transfer resistance vs. porosity is plotted in Fig. 10(b). 
The maximum performance is achieved when ¢ = 0.3. When the 
porosity is equal to zero, the total resistance drastically increases 
by almost ten times, since there are still some active TPB’s on 
the gas/cathode interface. On the other hand, the reaction sites 
completely disappear and the resistance goes to infinity when 
not any pore exists in the cathode. The above results have suffi- 
ciently demonstrated the porosity 0.5 is the upper limit and 0.3 
is the lower limit of porosity in a composite cathode. 
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Fig. 11. Effect of thickness on reaction rate and concentration overpotential: (a) 
local average reaction rate; (b) local average concentration overpotential. 


5.3.2. Layer thickness effect 

To discover the dependence of the reaction rate and overpo- 
tentials on the cathode thickness, several cases with different 
thicknesses are employed. Fig. 11 gives the comparison results. 
When the thickness is less than 30 um, the curve for the reaction 
rate is almost flat, i.e., almost all active bonds of the cathode are 
used; no matter how far they are from the gas/cathode interface. 
However, when the thickness is increased beyond 30 um, only 
the zone which is close to the electrolyte/cathode is well utilized 
for reaction as shown in Fig. 11(a). The same tendency is also 
demonstrated in Costamagna et al.’s study [4]. A parameter, i.e., 
“effectiveness factor à” is defined as a criterion to indicate what 
extent the electrode is utilized. When À is greater than 0.3 and 
one of two phases has much higher resistivity than the other 
phase, the major part of the electrode becomes useless and the 
activation overpotential (or reaction) is concentrated on the elec- 
trolyte/cathode interface. Moreover, with increasing thickness, 
the concentration loss will become more and more significant 
as shown in Fig. 11(b) due to longer diffusion path and poor 
gas supply. Comprehensively considering the above results, sug- 
gested that the cathode thickness should be less than 89 um here 


1000 
experimental data,y,=0.32 
800 i experimental data, y =0.4 
è continuum theory, y,=0.32 
600 3 --°-- continuum theory, y,=0.4 


simulated result, y =0.32 
simulated result, y =0.4 


400 


200 


Specific total resistance (m Q-cm’) 


0 10 20 30 40 50 60 70 
Thickness (um) 


Fig. 12. Comparison among the model prediction, experimental data and con- 
tinuum model. 


to achieve minimized concentration loss and uniform reaction 
rates. 


5.3.3. Combined effect of composition, particle size and 
thickness 

One of the important criteria indicating the electrode per- 
formance is the total specific resistance. Some experimental 
results have demonstrated the possibility of optimized geomet- 
rical parameters to minimize the resistance. In this section, 
we will consider the combined effects of composition, par- 
ticle size and thickness on this resistance. The validation of 
the model with the experimental data and other theoretical 
study is presented in Fig. 12. We have incorporated literature 
experimental parameters and some reasonable approximated 
values [4,29,30] for an EDB/Pt cathode into our model, e.g., 
€=0.4, T=900°C, y=60°, ¢=0.4, ry =0.1 um, Fio =0.1 wm, 
ooe1=2.5 x 10°Sm7!, 09,9 =20Sm7!. Kenjo et al. [29,30] 
does not give a relationship between overpotential and exchange 
current density such as Eq. (9), but an approximation for 
exchange current density, i.e., iọ =400 A m7? was reported by 
the same authors and this value has been adopted in our model. 
The predicted results are in good agreement with the experi- 
mental data by Kenjo et al. [29,30]. When the thickness exceeds 
a certain value L“, the specific polarization levels off, keep- 
ing an almost constant trend. This can be explained by that 
the electronic conductor is close to a pure conductor Oel > dio, 
the charge flows through an electronic conductor without any 
ohmic losses and the most of electrochemical reactions occurs 
in a narrow region next to the electrolyte/cathode with the 
thickness is increased as shown in Fig. 11(a). Therefore, the 
electrode resistance does not increase as the ohmic losses are 
negligible. Compared with the results from a continuum the- 
ory, there is a significant derivation as the electrode is very 
thin (<10 um). This derivation becomes more remarkable for 
the case of a lower electrode volume fraction g=0.32. The 
resistance calculated from our model is higher than that from 
Costamagna et al.’s model [4]. The continuum theory implies 
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Fig. 13. Effect of electrode fraction volume and particle radius on total resistance: (a) for different pe; (b) Ye] = 0.4 and different radius; (c) ge) = 0.2 and different 


radius. 


an assumption that only “percolation cluster” dominates the 
performance and the “isolated” and “dead-end” clusters are 
negligible. Therefore, for a very thin electrode in which that 
the percolating cluster almost vanishes, this theory tends to 
underestimate the resistance. For a thicker electrode, the con- 
sistency between the two models is satisfactory. For the case 
with e1 =0.32, our model reveals a smaller value of L”, while 
Costamagna et al.’s result indicates a larger L“ value and gradual 
decrease in resistance before Z<75 um. Simulated results suffi- 
ciently indicate that our model can well predict the experimental 
data. 

Fig. 13 presents the combined effects of electrode volume 
fraction, thickness and radius of particle on the total resis- 
tance. The continuum theory and the theory of coordination 
number demonstrate that the percolation threshold varies from 
0.294 to 0.709 as the radius ratio re1/fio equals to unity. Within 
these thresholds, reasonable performance of the cathode can be 
achieved. The resistance distribution within the cathode as a 
function of the electrode fraction is displayed in Fig. 13(a). The 
porosity, radius and electrode volume fraction are fixed at £ = 0.4, 
el = Fio = 1 wm and oe =0.5, respectively. When the electrode 
fraction is equal to 0.35, 0.4 and 0.5, respectively, the resis- 
tance maintains a constant after a certain thickness value and the 
cathode thickness corresponding to the optimum performance 
is in the range from 17.5 to 40 um. The maximal performance 
is achieved at pe =0.4, and Z=29.14 um. However, for cases 


with e1 < 0.3, the resistances become very large due to insuffi- 
cient reaction sites. If considering a cathode only composed of 
the LSM phase with parameters: Lz = 29.76 um, rel = Fio = 1 wm, 
€=0.4, the calculated total resistance is about 1695 mQ cm?. 
Compared this value with the data shown in Fig. 13(a), the 
advantage of a composite cathode is clear. From Fig. 13(b) and 
(c), it can be observed that with the decreasing of the radius, the 
optimal thickness corresponding to minimal total resistance also 
decreases when the electrode fraction is within the thresholds. 
This implies that the larger particle size, the thicker the cath- 
ode is needed to reach an optimum performance. While beyond 
the suggested percolation, i.e., ge] = 0.2, the optimum thickness 
does not exist for all three cases. Therefore, the performance of 
a cathode can be improved by decreasing the particle size and 
thus enhancement in the TPB length per unit volume. However, 
very small LSM or YSZ particles correspondingly increase the 
gas diffusion path and make the pore size become smaller, so the 
concentration overpotential will increase, especially for a much 
thicker cathode. 


6. Conclusions 


A 3D random micro-scale model taking into account of the 
gas diffusion and pore structure has been developed to eval- 
uate the mass transfer and electrochemical performance for a 
composite cathode in a solid oxide fuel cell. Then, detailed inves- 
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tigations are carried out to optimize the geometry parameters. 
A comparison with literature experimental data shows a good 
agreement. 

Model results demonstrate that the concentration loss is 
insignificant only when the porosity is in the middle range and 
the cathode is thin under a lower current density. The electro- 
chemical reaction behavior is strongly influenced by the porosity 
value. Only when the porosity is in the range from 0.25 to 0.45, 
the gas penetration depth is larger, reaction rate is relatively uni- 
form and almost all of the TPB sites become the active TPB 
sites. On the one hand, the larger thickness is favorable for the 
performance by increasing the number of the effective reaction 
sites and decreasing the specific resistance. While on the other 
hand, the thickness should be limited to less than a certain level 
in order to avoid excessive concentration loss and lower utiliza- 
tion. Smaller particle sizes reduce activation overpotential by 
increasing the TPB length per unit volume of the cathode. The 
minimum value in total specific resistance decreases as the parti- 
cle radius is decreased. But at the same time, too small particles 
will increase the gas diffusion resistance. In the present study, 
the best performance is achieved at Ye; = 0.4, and Z=29.14 um 
when £ =0.4. 

There is a good agreement between our predicted results 
with the continuum theory and experimental data. It shows that 
the predicted results on total specific resistance are reasonably 
larger than that obtained from the continuum theory for a very 
thin cathode. The presented model is still valid for the cathode 
with a lower ratio between the thickness and particle diameter 
compared with the continuum theory. Also, the current model 
removes some weaknesses in the Sunde’s model [18,19], such 
as ignoring the diffusion resistance and specific pore structure. 
Therefore, this model can be applied to much more general 
situations. 
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